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Abstract 

Background: The Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome is present in at least 1 out of 4,500 female live 
births and is the second most common cause for primary amenorrhea. It is characterized by vaginal and uterine 
aplasia in an XX individual with normal secondary characteristics. It has long been considered a sporadic anomaly, 
but familial clustering occurs. Several candidate genes have been studied although no single factor has yet been 
identified. Cases of discordant monozygotic twins suggest that the involvement of epigenetic factors is more likely. 

Methods: Differences in gene expression and methylation patterns of uterine tissue between eight MRKH patients 
and eight controls were identified using whole-genome microarray analyses. Results obtained by expression and 
methylation arrays were confirmed by qRT-PCR and pyrosequencing. 

Results: We delineated 293 differentially expressed and 194 differentially methylated genes of which nine overlap 
in both groups. These nine genes are mainly embryologically relevant for the development of the female genital 
tract. 

Conclusion: Our study used, for the first time, a combined whole-genome expression and methylation approach 
to reveal the etiology of the MRKH syndrome. The findings suggest that either deficient estrogen receptors or the 
ectopic expression of certain HOXA genes might lead to abnormal development of the female reproductive tract. 
In utero exposure to endocrine disruptors or abnormally high maternal hormone levels might cause ectopic 
expression or anterior transformation of HOXA genes. It is, however, also possible that different factors influence the 
anti-Mullerian hormone promoter activity during embryological development causing regression of the Mullerian 
ducts. Thus, our data stimulate new research directions to decipher the pathogenic basis of MRKH syndrome. 



Background 

The Mayer-Rokitansky-Kuster-Hauser (MRKH) syn- 
drome (OMIM 277000) is the second most common 
cause of primary amenorrhea and affects at least 1 in 
4,500 females. It is characterized by congenital absence 
of the uterus and the upper two thirds of the vagina in 
women with a normal female karyotype. As the ovaries 
are functional, women affected have physiological hor- 
mone levels and normal secondary sexual characteristics 
[1]. The MRKH syndrome may occur isolated (type I), 
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or can be associated with renal or skeletal malforma- 
tions, and, to a lesser extent, auditory and cardiac 
defects (type II) [2]. Although it is generally sporadic, 
familial clustering has been described, indicating a 
genetic cause [3]. Familial cases have been explained by 
autosomal dominant inheritance with incomplete pene- 
trance and variable expressivity or by small chromoso- 
mal aberrations undetectable in standard karyotypes [4]. 
However, the lack of families with informative genetic 
histories has not allowed the identification of any locus 
using standard genetic linkage analysis. Investigations 
have therefore used a candidate gene approach based on 
association with other genetic diseases or involvement 
during embryogenesis [4]. 
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The association of abnormalities in Mullerian duct 
(MD) development with renal, skeletal, cardiac and 
auditory defects suggests that crucial genes of fetal 
development and sex differentiation such as HOX, WNT 
and those encoding anti-Miillerian hormone (AMH) and 
its receptor are potential candidates [5,6]. The WNT 
genes control the production of a large family of pro- 
teins involved in intercellular signaling during embryo- 
genesis. Heterozygous mutations of WNT4 have been 
detected in a subgroup of patients, but these patients 
also show signs of hyperandrogenism [7-10]. HOX genes 
play key roles in body patterning and organogenesis, in 
particular during genital tract development and the dif- 
ferentiation of the kidneys and skeleton. Thus, expres- 
sion or function defects in one or several HOX genes 
may account for this syndrome. Furthermore, several 
hormones regulate physiological processes in the adult 
female reproductive tract by regulating HOX gene 
expression. Alterations in HOX gene expression that 
persist in the adult are a molecular mechanism by 
which endocrine disruptors may affect reproductive 
tract development [11,12]. However, structural abnorm- 
alities in HOX genes or in hormones regulating HOX 
expression have not been identified in women with 
MRKH syndrome until today [13-17]. As a third group 
of genes, AMH and its receptor have been regarded as 
causative factors in MRKH syndrome as AMH initiates 
MD regression in the 6th gestational week [18-22]. 
Mutation analyses of the AMH gene, however, did not 
support a link between MRKH syndrome and AMH yet 
[18,23]. Finally, mutations in other genes with a broad 
spectrum of activity during early development such as 
WT1, PAX2 and others have also been excluded in 
MRKH patients [4]. 

From the previous studies one can conclude that the 
targeted candidate gene approach has failed to decipher 
the causes of MRKH syndrome [2,5,24]. Recently, sev- 
eral recurrent copy number variants in patients with iso- 
lated and syndromic Mullerian aplasia have been 
described, but none of them was consistently found in a 
larger group of patients [25]. 

Cases of discordant monozygotic twins suggest that 
the involvement of epigenetic factors is more likely. Sev- 
eral studies have identified epigenetic differences, either 
for selected genes in monozygotic twins or in the overall 
epigenome [26,27]. 

We provide here the first study using a whole-genome 
approach to detect differences at the transcriptome and 
methylome level between MRKH patients and healthy 
controls. As integrated genomics becomes more and 
more important, the synergy between transcriptional 
and epigenetic gene regulation may be used to better 
understand the etiology of MRKH syndrome. 



Methods 

Patients 

This study was approved by the ethics committee of the 
Eberhard-Karls-University of Tuebingen. Between July 
2007 and December 2010, we had partly or completely 
excised 102 rudimentary uterine structures during 
laparoscopic-assisted neovagina in MRKH patients after 
informed consent was obtained [28]. As controls, we 
included 63 patients who underwent hysterectomy for 
benign disease in the same period. 

Microarray analysis was performed in eight patients 
and eight controls to detect differentially expressed 
genes and seven patients and seven controls from the 
same group to detect differentially methylated CpG 
sites. Of these patients, four had MRKH type I and four 
had MRKH type II, including three patients with skeletal 
malformations, and amongst these, one with Fallot's tet- 
ralogy and one with ureter abnormalities. None of the 
patients had MURCS association or other complex 
malformations. 

Analysis of serum samples at the time of surgery 
showed similar distribution between cycle phase one 
and two in the patients and control group. 

Tissue samples were examined histologically before 
RNA and DNA were isolated. All tissue samples in both 
groups consisted of more than 80% myometrium. 

RNA and DNA isolation 

The total RNA from myometrial pieces of rudimentary 
uterine tissue or normal uterus was isolated using the 
RNeasy® Mini Kit (Qiagen, Hilden, Germany). RNA 
quality was checked by a Lab-on-a-Chip-System Bioana- 
lyzer 2100 (Agilent, Boeblingen, Germany), and the con- 
centration was determined using a BioPhotometer 
(Eppendorf, Hamburg, Germany). DNA was isolated 
using the DNeasy® purification Kit (Qiagen, Hilden, 
Germany) according to protocol and the concentration 
was determined using a BioPhotometer (Eppendorf). 

Affymetrix microarray analysis 

Double-stranded cDNA was synthesized from 100 ng of 
total RNA and subsequently linearly amplified and bioti- 
nylated using the GeneChip® WT cDNA Synthesis and 
Amplification Kit (Affymetrix, Santa Clara, CA, USA) 
according to the manufacturer's instructions. 15 \ig of 
labeled and fragmented cDNA was hybridized to Gene- 
Chip® Human Gene 1.0 ST arrays (Affymetrix). Arrays 
were scanned using the GCS3000 Gene Chip scanner 
(Affymetrix) and AGCC 3.0 software. Scanned images 
were inspected visually to check for hybridization arti- 
facts and proper grid alignment and analyzed with 
Expression Console 1.0 (Affymetrix) to generate report 
files for quality control. 
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Quantitative real-time PCR 

Relative expression of selected mRNA targets was deter- 
mined by quantitative real-time PCR (qRT-PCR). 250- 
500 ng of total RNA was reverse transcribed using a 
QuantiTect Reverse Transcription Kit (Qiagen) accord- 
ing to the manufacturer. cDNA was diluted 1:10 before 
PCR amplification or preamplified using the TaqMan 
PreAmp Master Mix Kit (Applied Biosystems, Carlsbad, 
CA, USA) according to the manufacturer's protocol and 
diluted 1:20 for the subsequent PCR analysis. Primers 
were designed with Primer3 or PrimerBlast (http://bio- 
tools.umassmed.edu/bioapps/primer3_www.cgi, http:// 
www.ncbi.nlm.nih.gov/tools/primer-blast/) and synthe- 
sized by Metabion (Metabion, Martinsried, Germany). 
Table 1 gives a list of PCR targets and primers. 

Real-time detection of specific PCR products was per- 
formed on a LightCycler480 (Roche, Penzberg, Germany) 
with 5 uL of 2x QuantiTect SYBR Green PCR Kit (Qia- 
gen). The PCR reaction was initiated by a 10 min hot start 
followed by 45 cycles of 95°C for 20 s, 58°C for 40 s, and 
72°C for 20 s. Each PCR reaction was performed in three 
technical replicates. PCR efficiency was calculated from 4- 
or 5 -fold serial dilutions of an equal mixture of all cDNAs 
using the following equation: E = 10 A [-l /slope] [29]. To 
calculate the relative expression of each target, the raw Cp 
values were imported into qBase [30] . Three suitable refer- 
ence genes (PDH, SDHA, PGRMC1 or HISPPD1) were 
selected according to their M-values and used for normali- 
zation of the qRT-PCR reactions [31,32]. 

Illumina methylation array analysis 

200-500 ng DNA was bisulfite-converted using the EZ 
DNA Methylation Kit (Zymo Research, Orange, CA, 
USA) according to the manufacturer's protocol. 



Seven patient and seven control probes were evaluated 
for genome-wide promoter methylation using the Illu- 
mina Infinium HumanMethylation27 BeadArray. 

After bisulfite conversion, each sample was whole-gen- 
ome amplified (WGA) and enzymatically fragmented. 
The bisulfite-converted WGA-DNA samples were puri- 
fied and applied to the BeadChips. Allele-specific primer 
annealing is followed by single-base extension using 
DNP- and Biotin-labeled ddNTPs. DNA methylation 
values, described as beta values, are recorded for each 
locus in each sample via BeadStudio software [33]. 

Differential methylation was assessed by subtracting 
the mean methylation level (beta value) of the patient 
group from the mean beta value of the reference group 
using BeadStudio software. 

Pyrosequencing methylation analysis 

PCR and sequencing primer were designed by Qiagen 
using the PyroMark Assay Design Software 2.0 and are 
shown in Table 2. For each gene, we selected the 
sequence of the CpG island region which had previously 
been identified as differentially methylated in the array 
experiments. One to five adjacent CpG sites were ana- 
lyzed for each CpG island. Five patient and five control 
samples were included in the experiment. 200 ng of iso- 
lated DNA was bisulfite-converted using the EZ DNA 
Methylation Kit (Zymo). 24 ng of bisulfite-treated DNA 
were amplified in reaction mixture containing forward 
and reverse primer, 1 U of HotStarTaq DNA Polymer- 
ase (Qiagen), 200 umol each of dNTP/1 and nuclease- 
free water. The same cycling conditions were used for 
all assays: denaturing at 95°C for 15 min; 35 cycles at 
95°C for 10 s, at 59, 3°C for 30 s, and at 72°C for 30 s; 
an additional elongation step was performed at 72°C for 



Table 1 qRT-PCR targets with corresponding Affymetrix probeset ID and primer used for amplification 



Target Affymetrix cluster ID Forward primer Reverse primer 



H0XA5 


8138735 


CGCCCAACCCCAGATCTA 


GGCCGCCTATGTOTCATG 


H0XA9 


8138749 


GOTGTGGTOTCCTCCAGT 


CCAGGGTCTGG I G I I I IGIA 


PGR 


7951165 


TGGTGmGGTCTAGGATGGA 


GGATCTGCCACATGGTAAGG 


ESR1 


8122840 


GCAGGGAGAGGAGmGTGT 


CAGGACTCGGTGGATATGG 


OXTR 


8085138 


GCACGGTCAAGATGACmC 


GCATGTAGATCCAGGGGTO 


PEG 10 


8134339 


GACCCCATCOTCCTGTOT 


GCTOACTOTGTGGGGATG 


MFAP5 


7960919 


TGCTCTCGTOTGTCTGTAAGG 


ACAGGGAGGAAGTCGGAAGT 


IRS1 


8059470 


GmCCAGAAGCAGCCAGAG 


GGAAGATATGAGGTCCTAGTOTGA 


IRS2 


7972745 


OTOTGTCCCACCAOTGA 


CAGTGCTGAGCG I C I IU I I I 


IGF2 


7937772 


ACACCCTCCAGTOGTCTGT 


CGGAAACAGCACTCCTCAA 


WISP2 


8062864 


GCGACCAACTCCACGTCT 


GTCTCCCCTOCCGATACA 


CDH5 


7996264 


ACAACGAGGGCATCATCAA 


AATGACCTGGGCTCTGmC 


SDHA 


8104166 


AGAAGCCCmGAGGAGCA 


CGATOCGGGTCTATATOCAGA 


PDHB 


8088384 


GAGGCTGGCCACAGmG 


GAAATOAACGCAGGACC^ 


PGRMC1 


8169617 


GGTGTOGATGTGACCAAAG 


TGAGGTCAGAAAGGTCATCGT 


HISPPD1 


8169617 


TCCATCATCTGACGTOCAC 


TGGTGTOGGAGGATCmG 
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Table 2 Pyrosequencing targets with corresponding lllumina probeset ID and primer used for amplification 



Target lllumina Target ID Forward primer 



Reverse primer 



Sequencing primer 



WISP2 cg03562120 

H0XA5 cg02248486 

H0XA9 cg27009703 

WT1 cg25094569 

GATA4 eg 17795240 



GTGTGTGmGGGAGTGAm 

GGAAmTGA^mATAATOT 
GTAATOGTAG^ 

Bi o-GTGGTG ATGGTGGTGGTATAT 



Bio-CTCATATCCCCTACAAAACCAACmAA 
Bio-AACCACAAATCAAACACACATATCA 



AOTCAACCCCTACAAOTCCAATCCA 
Bio-TGGATGTGA^GGGATAGGT CCCA^mAAAACCAAACCAmAACT 



AAGGATOGmAGGGAGAGmG^G Bio-TAAAAmCACCATATOACCAAAAACT 

CCTAACC™ 



GmGGGAGTGA^ATAGTOT 

AAmTGA^mATAATOTGTAA 
TOGTAG 

TCAACCCCTACAAOTCCAATCCA 
AAA 

A^mAAAAAATAAACAACCTO 
TCTATC 

GGmAGGGAGAGmG^G 



3 min. Gel electrophoresis was carried out on all PCR 
products. All PCR reactions included a no-template con- 
trol and four standardized methylation controls (0%, 
30%, 70%, and 100% methylated DNA). Pyrosequencing 
was carried out using the Pyrosequencer PSQ 96 MA 
(Biotage AB, Uppsala, Sweden). Results were automati- 
cally analyzed using the PSQ 96MA 1.0 software (Allele 
Quantification mode). 

Pathway analysis 

Gene regulation networks were generated using Ingenu- 
ity Pathway analysis software (ingenuity® systems, 
http://www.ingenuity.com). The dataset with differen- 
tially regulated transcripts and their corresponding 
expression and methylation values were uploaded into 
the application. The genes were overlaid onto a global 
molecular network developed from information in the 
Ingenuity Pathways Knowledge Base. Networks of these 
focus genes were then algorithmically generated based 
on their connectivity. All edges are supported by at least 
one reference from the literature, from a textbook, or 
from canonical information stored in the Ingenuity 
Pathways Knowledge Base. Genes from the microarray 
dataset that met the fold change cutoff of 1.5 and that 
were associated with a relevant pathway in the Ingenuity 
Pathways Knowledge Base were included in the analysis. 

Statistical methods 

With quantitative RT-PCR normalized values were 
obtained from qBase. The ratio of means of patients 
divided by means of controls was then calculated and 
the log2 of this ratio is shown. For the Affymetrix 
microarray analysis the means of patients were divided 
by the means of controls and the log2 was shown as log 
fold change. Measurement errors were calculated using 
Gaussian error propagation. As usual it is assumed that 
normalized values obtained from qBase are lognormal 
distributed and a t-test (i.e. the Welch test assuming 
unequal variances) was applied for each gene to investi- 
gate the difference between patients and controls in 
qRT-PCR. The same was done for expression array data. 
A significance level of 5% was chosen. The percentage 



of methylated cytosines was obtained from pyrosequen- 
cing methylation analysis. The difference means of 
patients minus means of controls are shown. The same 
was done for the average beta values. 

Results 

Microarray expression and methylation analysis 

To identify changes in the expression level of putative 
candidate genes, we performed microarray analysis with 
Affymetrix Human Gene 1.0. Analysis using the 
ArrayAssist 4.0 software identified 293 transcripts differ- 
entially expressed between tissue samples of MRKH 
patients and controls. Of these transcripts, 161 were 
upregulated and 132 downregulated with a fold change 
of at least 1.5 and a p-value of less than 0.05 (Table 3). 

Pathway analysis revealed genes relevant in the embry- 
ological development of the genital tract, including 
HOXA genes and hormone receptors. 

The delineation of regional DNA methylation patterns 
has important implications for understanding why cer- 
tain regions of the genome can be expressed in specific 
developmental contexts and how epigenetic changes 
might enable aberrant expression patterns and disease 
[27]. We therefore decided to compare whole-genome 
expression and methylation patterns in uterine rudi- 
ments of MRKH patients compared to control uteri. To 
achieve this, we performed lllumina HumanMethyla- 
tion27 BeadArrays and overlaid both datasets. The ana- 
lysis using the BeadStudio software identified 194 



Table 3 Differential expression of genes and methylation 
of CpG-sites in MRKH patients compared to controls in 
numbers 



Difference in patients from controls 


n 


Differential expression, total 


293 


Downregulated 


132 


Upregulated 


161 


Differential methylation, total 


194 


Hypomethylated 


116 


Hypermethylated 


78 


Overlap 


9 
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differentially methylated CpG sites in specific CpG 
islands. Of these sites, 78 were hypermethylated and 116 
hypomethylated (Table 3). 

Nine genes were detected in both datasets (HOXA5, 
HOXA9, WISP2, CDH5, PEG 10, MFAPS, LRRC32, 
RALGPS2, SMPD3); these are termed 'overlap genes' 
(Table 4). CpG sites within these genes were either 
hypermethylated and the genes underexpressed or hypo- 
methylated and overexpressed, except for one gene. 

Of the nine overlap genes, six (CDH5, MFAPS, WISP2, 
HOXA5, PEG 10, HOXA9) were included in the subse- 
quent analyses and experiments as they are known to be 
relevant to the embryological development of the female 
genital tract. 

Network and pathway analysis 

Ingenuity Pathways analysis software (Ingenuity Sys- 
tems) was used to examine the connection between the 
differentially methylated CpG sites and differentially 
expressed genes. As shown in Figures 1 and 2, differen- 
tially expressed genes and differentially methylated sites 
can be assigned to basic functions relevant to cell and 
tissue development and proliferation, cell-to-cell signal- 
ing and interaction, cellular movement, cancer, endo- 
crine and reproductive disorders, and others (Figure 1 
and 2). 

Figure 3 shows a network of differentially regulated 
overlap genes and other relevant genes. The network 
was created by fusing the expression and methylation 
datasets from the microarray experiments. This specific 
network was selected because of the known relevance of 
the genes included during embryological development 
and during functional changes of the female reproduc- 
tive tract. Interactions between embryologically relevant 
genes, including HOXA genes and hormone receptors, 



can be clearly detected. Interacting genes are either dif- 
ferentially expressed, carry differentially methylated CpG 
sites, or both. The gene regulation network contains 15 
differentially regulated genes, seven downregulated 
(green icons), two of these also with hypermethylated 
CpG sites (dark grey), and seven upregulated (red 
icons), three of these also with hypomethylated CpG 
sites (light grey). Four genes contained differentially 
methylated CpG sites without being differentially 
expressed (three hypo- and one hypermethylated). Eight 
genes are supplemented in the network to complete the 
interactions. Eleven of the transcripts shown were used 
for qRT-PCR validation and five for pyrosequencing. 

Figure 4 shows the relation between expression and 
methylation in overlap genes later selected for valida- 
tion. CpG sites in genes were either hypomethylated 
with genes overexpressed or hypermethylated with 
genes underexpressed, except for CDH5. 

Validation of expression differences by qRT-PCR 

We chose nine key players in the interaction networks 
for independent verification by qRT-PCR. The genes 
were selected because of their known relevance during 
embryological development. Three suitable reference 
genes (PDH, SDHA, PGRMC1 or HISPPD1) were 
selected according to their M-values and used for nor- 
malization of the qRT-PCR reactions. We were able to 
validate all nine genes. The results of the qRT-PCR (Fig- 
ure 5) showed 100% validation efficiency in comparison 
to the expression data of the microarray experiment, 
although statistical significance was not always found. 

Validation of methylation differences by pyrosequencing 

Five embryologically important genes were chosen for 
validation of the methylation array experiments by 



Table 4 Overlap genes: names of differentially expressed genes that contain differentially methylated CpG-sites 



Probe set ID, 
human 1.0 
genechip array 


Probe set ID, 

27human 

methylation 


Gene title 


Gene 
symbol 


Fold change 
human 1.0 
genechip array 


P- 

value 


Diff methyl 
27human 
methylation array 


P- 

value 


8138735 


cg02248486 


homeobox A5 


HOXA5 


1.9 


0.00036 


-0.33 


0.00015 


8138749 


cg26521404 


homeobox A9 


HOXA9 


1.5 


na 


-0.23 


0.00047 


8062864 


cg03562120 


WNT1 -inducible signaling pathway 
protein 2 


WISP2 


-1.7 


0.00148 


0.14 


0.00096 


7996264 


cg22319147 


cadherin 5, type 2, VE-cadherin 
(vascular epithelium) 


CDH5 


-1.6 


0.00503 


-0.12 


0.00598 


8134339 


eg 191 07595 


paternally expressed 10 


PEG10 


1.8 


0.00290 


-0.12 


0.01596 


7960919 


cg1 581 5843 


Microfibrillar-associated protein 5 


MFAP5 


-2.1 


0.04635 


0.16 


0.02451 


7950555 


cg20899321 


Leucine-rich repeat containing 32 


LRRC32 


-1.6 


0.00435 


0.12 


0.00045 


7907657 


eg 10559803 


Ral GEF with PH domain and SH3 
binding motif 2 


RALGPS2 


2.0 


0.00385 


-0.11 


0.00345 


8002249 


eg 172 17677 


sphingomyelin phosphodiesterase 3, 
neutral membrane (neutral 
sphingomyelinase II) 


SMPD3 


1.6 


0.00001 


-0.10 


0.01602 
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Figure 1 Assignment of differentially expressed genes to functional groups. Assignment of all differentially expressed genes to functional 
groups. The analysis was done with Ingenuity Pathways analysis software. Fischer's exact test was used to test for significance (shown as bars), 
determining the probability that each biological function assigned to the network is due to chance alone. 



pyrosequencing. We selected one to five CpG sites 
within one specific CpG island per gene for analysis 
(Figure 6). The CpG islands were selected according to 
the differential methylation in the preceding array 
experiments. The differential methylation status of the 
array experiments was confirmed for all five CpG 
islands within the WISP2, HOXAS, HOXA9, GATA4 



and WT1 genes, thus validation efficiency was again 
100% (Figure 7). 

Discussion 

Although MRKH syndrome is a congenital disorder 
most patients are not diagnosed until puberty. Using a 
candidate gene approach, the underlying cause has so 




llllllllliii 



ment 


ion 


on 


ion 


functi 




o 


1 


ifer 


tera 


C 


arol 


g 


Cellular d 


13 






evelopme 


growth ai 


signaling 




3 


lar 






~s. 


r>llll 






>. 

00 


U 






cal 




U 




f 







C -3 



C3 

0 



Figure 2 Assignment of all genes containing differentially methylated CpG-sites to functional groups. Fischer's exact test was used to 
test for significance (shown as bars), determining the probability that each biological function assigned to the network is due to chance alone. 
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Figure 3 Network of differentially expressed genes and genes containing differentially methylated CpG-sites. The fold change of 
regulation in the microarray analysis, (the p-value of the significance analysis) and the percentage of differential methylation are listed below the 
symbols. For the purposes of simplification, only selected known gene-to-gene interactions are shown. Continuous arrow lines show direct 
interactions between genes and broken arrow lines show indirect interactions. 



CDH5 MFAP5 WISP2 H0XA5 PEG10 H0XA9 



cu 



c 
o 

GO 
>, 
-C 

CD 

E 
c 

CD 
O 

c 

CD 

I 



CO 

d 



CM 

d 



o 
d 



cm 
d 



CO 

d 



d 



LO 

d 



LO 

d 



TO 



Q_ 

X 
CD 
CD 
CD 
C 
CD 
-C 
O 

o 

4— 

CD 
O 



Figure 4 Relation between expression and methylation in overlap genes. Methylation (left) and expression (right) of overlap genes in 
MRKH patients compared to control group. Methylation is shown as difference of average beta of patients minus controls. Expression differences 
are shown as log fold change of patients divided by controls. Bars indicate the measurement error. 
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Figure 5 Validation of differential expression by qRT-PCR. Expression of overlap genes in MRKH patients compared to control group 
measured by array (left) and qRT-PCR (right). Array data are shown as log fold change of patients divided by controls, qRT-PCR data are shown 
as log ratio of patients divided by controls. Bars indicate the measurement error and stars indicate a significant difference between patients and 
controls (t-test, p-values <0.05). 



far not been identified [5]. Recently, a high incidence of 
recurrent copy number variants in patients with isolated 
and syndromic Mullerian aplasia has been described, 
but none of them was consistently found in a larger 
group of patients [25]. Cases of discordant monozygotic 



twins suggest that the involvement of epigenetic factors 
is more likely. The present study was the first to use a 
whole-genome approach to identify relevant genes, 
including differential expression and methylation. This 
allowed us to create a complex network of genes, which 
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Figure 7 Validation of differential methylation by pyrosequencing. Methylation array (left) and pyrosequencing (right) of overlap genes in 
MRKH patients compared to control group. Both are shown as difference of patients minus controls. Bars indicate the measurement error. 



provided insight into possible mechanisms underlying 
MRKH syndrome. Our data indicate that different 
potential mechanisms are possible. 

Deficiency of hormone receptors 

The overexpression of ESR1 and PGR in rudimentary 
uterine tissue from MRKH patients could be explained 
by a deficiency of these hormone receptors. The local 
overexpression may be the result of a positive feedback 
mechanism well known from other hormonal regulatory 
loops. Similar hypotheses have been postulated before 
but firm scientific evidence has not yet been obtained. 
This is the first study with an experimental setting that 
supports the hypothesis of locally deficient hormone 
receptors in MRKH syndrome. 

As early as 1910, Kuster explained the MRKH syn- 
drome by regression of the Miillerian duct (MD) [34]. 
According to Ludwig, the MRKH syndrome results from 
non-fusion of the MD with the Wolffian duct (WD) 
[35,36]. Because the embryo is under the influence of 
maternal hormones, he suggested that both the non- 
fusion of the MD with the WD and rudimentary devel- 
opment of the vagina are caused by a deficiency of 
gestagen or estrogen receptors [35,36]. Estrogens are 
necessary for the embryonic development of the female 
reproductive tract. The special role of ESR1 in female 
reproductive tract development has been demonstrated 
by disrupting the corresponding gene in the mouse, 
resulting in hypoplastic uterine and vaginal tissue [37]. 

Influence of estrogen on AMH and its receptor 

One of the first hypotheses for the underlying cause of 
MRKH syndrome was an activating mutation of either 
the gene for the AMH receptor (AMHR), resulting in 



the inappropriate production of AMH, or the receptor 
itself. Schmid-Tannwald and Hauser proposed the 
regression of the MD due to a temporary secretion of 
AMH during the first fetal weeks [38]. Depending upon 
the amount secreted, a greater or lesser portion of the 
MD regressed. Nevertheless, mutation analyses of the 
AMH gene did not support a link between the MRKH 
syndrome and AMH at the genome level. Also, AMH 
protein levels in plasma and peritoneal fluid from 
MRKH patients were equivalent to control individuals 
[18,23]. Our study confirmed the data given, as we did 
not see any persistent differential expression or methyla- 
tion patterns in the AMH or AMHR genes in adolescent 
MRKH patients. However, low or baseline AMH levels 
in a female adolescent may not necessarily be correlated 
with the patient's early embryonic exposure to AMH 
signaling [23]. 

It has been reported that estrogen regulates AMH 
expression [39]. The constant overexpression of ESR1 
found in the rudimentary uterine tissue in our patients 
or the in utero exposure to abnormally high maternal 
levels of E 2 could lead to an increasing AMH promoter 
activity during embryological development of the female 
genital tract causing uterine and vaginal aplasia. 

Transcription factors involved in primary sex determi- 
nation are also recruited as important regulators of 
AMH transcription. A common regulatory factor impor- 
tant for transcription of AMH genes is WT-1. WT-1 is 
essential for the embryonic development of the kidneys 
and gonads. GATA4 appears to play a predominant role 
in sex determination and sex differentiation via AMH 
gene regulation [40]. In our study CpG sites within 
WT1 and GATA4 were both hypomethylated compared 
to the control tissue. This could be a sign of stable 
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activation leading, at least during embryological develop- 
ment, to activation of the AMH gene and thus partial 
regression of the MD. 

Influence of endocrine disruptors on the Mullerian duct 
and on the expression of HOX genes 

Chemical compounds homologous to steroids can act as 
agonists or antagonists in fetuses exposed to them [37]. 
The involvement of ED with estrogen-like functions 
would be an explanation for the findings in MRKH syn- 
drome, although the analyses of pregnancy histories 
available in all our cases so far failed to identify any 
clear association with drug use, illness or exposure to 
known substances. 

Several examples of a negative impact of estrogens on 
uterine development are known. Transient exposure of 
the neonatal ewe to estrogens during critical periods 
specifically perturbs normal development of the uterus 
[41]. Estrogen inhibited caudal progression of develop- 
ing MD in the turtle. It has been shown to block devel- 
opment of the MD when applied before the start of 
differentiation, and the length of the MD varied with 
the time point treatment was given [42]. 

Epithelial-mesenchymal differentiation in the murine 
MD is regulated by WNT signaling correlated with 
expression of HOX genes. Several nuclear hormonal 
receptors regulate the expression of multiple HOX 
genes. When a HOX gene is mutated, the body segment 
where it is normally expressed typically develops charac- 
teristics of the segment anterior to it, an effect known 
as anterior transformation. In contrast to Drosophila, in 
vertebrates, targeted mutation in a single HOX gene 
usually causes only a subtle transformation. This is 
because of genetic duplication and functional redun- 
dancy of adjacent genes [11]. The hypomethylation of 
specific CpG sites and corresponding overexpression of 
HOXA9 could be due to either exposure to a substance 
similar to diethylstilbestrol (DES) in utero or a deficient 
HOXA10 causing anterior transformation. 

It is known that HOXA9 is expressed at high levels in 
areas destined to become the fallopian tube, HOXA10 is 
expressed in the developing uterus, HOXA11 is 
expressed in the primordia of the lower uterine segment 
and cervix, and that HOXA13 is expressed in the ecto- 
cervix and upper vagina. This expression pattern has 
been preserved in mice and humans [43]. Microarray 
analysis has shown organ-specific changes in gene 
expression profiles in the oviduct, uterus, and vagina 
after DES exposure. Changes in HOX and WNT expres- 
sion might lead to abnormalities of segment-related 
positional identity in the upper part of the MD after 
DES exposure [43,44]. Sex steroids have been investi- 
gated in the regulation of the HOX genes at the 5end of 
the cluster that determine posterior development, 



including development of the reproductive tract [45-48]. 
Both, HOXA10 and HOXA11 expression, is upregulated 
by 17E-estradiol and progesterone. Changes in HOX 
gene expression are a potential marker for the effects of 
in utero drug use that may become apparent only at late 
stages of development [12]. In utero, DES exposure 
shifts HOXA9 expression from the oviducts to the 
uterus and decreases HOXA10 as well as HOXA11 
expression of the uterus causing a T-shaped' uterus 
with a tube-like phenotype [47]. In human uterine and 
cervical cell cultures, DES has induced HOXA9 or 
HOXA10 gene expression [12]. 

Continued HOX gene expression in the adult has been 
described in the reproductive tract and may be a 
mechanism to retain developmental plasticity [49]. Spe- 
cifically HOXA10 and 11 are expressed in the endome- 
trium and their expression varies in a menstrual cycle- 
dependent manner. Although no women with mutations 
in HOXA10 and HOXA11 have been described, patients 
with lower implantation rates have lower HOXA10 and 
HOXA11 expression in the secretory phase [43]. 

In addition to HOXA9, specific CpG sites in HOXA5 
were hypomethylated and the gene overexpressed. This 
gene has a crucial role in the specification of the cervical 
and upper thoracic region of the skeleton. Its correct 
expression is important for the proper patterning of the 
embryo. It has been shown that ectopic HOXA5 expres- 
sion results in abnormal differentiation [50]. In a similar 
manner, ectopic HOXA5 expression at the 5'end of the 
cluster might prevent normal differentiation of the MD 
or even regression. HOXAS is known as a transcriptional 
regulator of multiple target genes, two of which are p53 
and the progesterone receptor (PGR). The overexpression 
of PGR in patients may be induced directly by the over- 
expression of HOXAS [51]. 

Finally, neonatal DES exposure is also known to cause 
overexpression of IRS-1 and IGF2, both of which are 
included in our network [52]. 

Impact of WNT genes on uterine development 

The WNT genes and products form the WNT signaling 
pathway which controls developmental processes. Only 
the WNT4 gene has been clearly implicated in atypical 
MRKH syndrome before [9]. The phenotype of WNT9b 
mutants can be rescued by activation of WNT1 in the 
WD, identifying the canonical WNT pathway as a deter- 
minant signaling process in MD elongation [37]. A 
recent study excluded mutations in the coding 
sequences of WNT4, WNT5A, WNT7A and WNT9B in 
11 MRKH patients [6]. In our study, CpG sites within 
the WISP2 (WNT1 inducible signaling pathway protein 
2) gene were hypermethylated and the gene underex- 
pressed in rudimentary uterine tissue, thus the relevance 
is not clear yet. WISP2 (CCN5), a gene that is important 
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in smooth muscle cell proliferation and migration, is an 
estrogen-induced gene in the uterus [53]. 

Conclusion 

We were able to draw important conclusions from our 
study, the first to compare rudimentary uterine tissue 
from MRKH patients and uterine tissue from healthy con- 
trols. GATA4, WT1 and constant overexpression oiESRl 
might increase AMH promoter activity during embryolo- 
gical development, resulting in partial regression of the 
MD. Involvement of endocrine disruptors (ED) with estro- 
gen-like functions might mimic this effect. The deficiency 
of hormone receptors may result in their overexpression 
and cause both the non-fusion of the MD with the WD 
and rudimentary development of the vagina. 

The hypomethylation of specific CpG sites and the 
corresponding overexpression of HOXA9 may be due to 
either exposure to a substance similar to DES in utero 
or deficient HOXA10 causing anterior transformation. 
Ectopic HOXAS expression at the 5'end of the cluster 
might prevent normal differentiation of the MD. 

Using the synergetic approach of transcriptional and 
epigenetic regulation, our study has, for the first time, 
provided a deeper insight into the etiology of congenital 
vaginal and uterine aplasia, and has significantly 
advanced the explanation of MRKH syndrome. Further 
investigations will show which of our hypotheses is cor- 
rect, but it is already clear that hormone receptors and 
HOX genes appear to play a major role and should be 
in focus of further examinations. 
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